//
// Created by hitfyd on 18-4-26.
//

/**
 * 项目需要使用GLPK线性规划开发库，GLPK使用可以参考官方文档http://www.gnu.org/software/glpk/，下载源代码后查看doc文件夹中的说明
 *
 * 环境配置：
 * Debian/Ubuntu终端输入：
    sudo apt install libglpk-dev glpk-utils
 *
 * 编译配置：
 * 在CMakeLists.txt中添加：
    add_executable(xxx xxx.cpp)     #xxx.cpp即为使用了GLPK的源文件，xxx为可执行文件名称
    link_libraries(glpk)
    target_link_libraries( xxx glpk )       #xxx为可执行文件名称
 */
#include <bits/stdc++.h>
#include "glpk.h"

using namespace std;

const int XNum = 1000, FNum = XNum;
vector<pair<int, int>> X;//有限集X，每个元素代表一个点并按照在向量中的顺序编号代表
vector<set<int>> F, C;//子集合集F，每个子集中存储点的编号;子集合集C，即为最后的输出
set<int> selectset, unselectset;//已选择点的集合和尚未被选择点的集合
const int subsetsize = 20;//子集合包含点的最大个数

//初始化已选择点的集合（清空）和尚未被选择点的集合（全集，等于X）
void initSelected(){
    selectset.clear();
    unselectset.clear();
    for (int i = 0; i < XNum; i++)
        unselectset.insert(i);
}

void createRandomNodes() {
    set<pair<int, int>> nodeset;
    int times = 0;
    while (nodeset.size() < XNum) {
        nodeset.insert(pair<int, int>(rand() % 100, rand() % 100));
        times++;
    }
    X.clear();
    copy(nodeset.begin(), nodeset.end(), back_inserter(X));
    printf("The times of generating random points:%d\n", times);

    F.clear();
    set<int> subset;
    initSelected();
    int n, x;
    //选出S0中的20个元素
    while (subset.size() < subsetsize) {
        int i = rand() % XNum;
        subset.insert(i);
        selectset.insert(i);
        unselectset.erase(i);
    }
    F.push_back(subset);
    //找到集合
    while (unselectset.size() >= subsetsize) {
        n = rand() % subsetsize;
        if(n == 0)
            x = 0;
        else
            x = rand() % n;
        subset.clear();
        while (subset.size() < x) {
            int i = rand() % XNum;
            if (unselectset.find(i) != unselectset.end()) {
                subset.insert(i);
                selectset.insert(i);
                unselectset.erase(i);
            }
        }
        while (subset.size() < n) {
            int i = rand() % XNum;
            if (selectset.find(i) != selectset.end())
                subset.insert(i);
        }
        F.push_back(subset);
    }
    F.push_back(unselectset);
    while (F.size() < FNum) {
        subset.clear();
        n = rand() % subsetsize;
        while (subset.size() < n) {
            subset.insert(rand() % XNum);
        }
        F.push_back(subset);
    }
}

void printSet(set<int> &s) {
    printf("set size:\t%ld\t", s.size());
    set<int>::iterator it = s.begin();
    while (it != s.end()) {
        printf("%d, ", *it);
        it++;
    }
    printf("\n");
}

void printVectorSet(vector<set<int>> &F) {
    printf("vector size:\t%ld\n", F.size());
    vector<set<int>>::iterator it = F.begin();
    while (it != F.end()) {
        printSet(*it);
        it++;
    }
}

bool verify(){
    set<int> allset;
    for(int i = 0; i < C.size(); i++){
        set<int>::iterator it = C[i].begin();
        while (it != C[i].end()) {
            allset.insert(*it);
            it++;
        }
    }
    if(allset.size() == XNum)
        return true;
    return false;
}

//权值为1的加权集合覆盖问题的线性规划
void Set_Cover_Approx(){

    glp_prob *lp;
    int Size = XNum * XNum;
    int *ia, *ja;
    ia = (int *)calloc(1 + Size, sizeof(int));
    ja = (int *)calloc(1 + Size, sizeof(int));
    double *ar, z, x[1 + XNum];
    ar = (double *)calloc(1 + Size, sizeof(double));
    lp = glp_create_prob();
    glp_set_prob_name(lp, "sample");
    glp_set_obj_dir(lp, GLP_MIN);
    glp_add_rows(lp, XNum);
    for(int i = 1; i <= XNum; i++){
        glp_set_row_name(lp, i, to_string(i).data());
        glp_set_row_bnds(lp, i, GLP_LO, 1.0, 0.0);
    }
    glp_add_cols(lp, XNum);
    for(int i = 1; i <= XNum; i++) {
        glp_set_col_name(lp, i, "x1");
        glp_set_col_bnds(lp, i, GLP_LO, 0.0, 0.0);
        glp_set_obj_coef(lp, i, 1.0);
    }
    int ne = 0;
    vector<int> XF;
    for(int i = 0; i < XNum; i++){
        int num = 0;
        for(int j = 0; j < F.size(); j++){
            if(F[j].find(i) != F[j].end()) {
//                printf("%d, %d\n", i, j);
                num++;
                ne++;
                ia[ne] = i + 1, ja[ne] = j + 1, ar[ne] = 1.0;
            }
        }
        XF.push_back(num);
    }
    glp_load_matrix(lp, ne, ia, ja, ar);
    glp_simplex(lp, NULL);
    z = glp_get_obj_val(lp);
    for(int i = 1; i <= XNum; i++)
        x[i] = glp_get_col_prim(lp, i);
    glp_delete_prob(lp);
    printf("z = %g\n", z);
    for(int i = 1; i <= XNum; i++)
        printf("x%d = %g\t", i, x[i]);

    C.clear();
    for(int i = 1; i <= F.size(); i++){
        if(x[i] >= (double)1 / XF[i-1])
            C.push_back(F[i-1]);
    }
}

int main() {
    clock_t startTime, endTime;
    createRandomNodes();
    startTime = clock();
    Set_Cover_Approx();
    endTime = clock();
    printVectorSet(F);
    printf("\n\n\n\n\n");
    printVectorSet(C);
    printf("%ld\t%ld\t%ld\n", startTime, endTime, CLOCKS_PER_SEC);
    printf("Algorithm run time:\t%lfms\n", (double)(endTime - startTime) / CLOCKS_PER_SEC * 1000);
    printf("C集合大小： %ld,\tVerify Result: %s\n", C.size(), verify() == true ? "True" : "False");
    return 0;
}